Covariant Wigner Function Approach to the Boltzmann 
Equation for Mixing Fermions in Curved Spacetime 



Kazuhiro Yamamoto 
Department of Physics, Hiroshima University, 
Higashi-Hiroshima 739-8526, Japan 

Based on the covariant Wigner function approach we derive the quantum Boltzmann equation 
for fermions with flavor mixing in general curved spacetime. This work gives a rigorous theoretical 
framework to investigate the flavor oscillation phenomena taking the gravitational effect into ac- 
count. It is shown that the Boltzmann equation of the lowest order of the expansion with respect to 
h reproduces the previous result which was derived in the relativistic limit on the Minkowski back- 
■ ground spacetime. It is demonstrated that the familiar formula for the vacuum neutrino oscillation 

f"^ | can be obtained by solving the Boltzmann equation. Higher order effects of the ^-expansion is also 

f^) , briefly discussed. 
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The Boltzmann equation, which describes non-equilibrium systems of many particles, is one of the most fundamen- 
tal to investigate physical processes in the early universe. For example, the statistics of the angular anisotropics in 
the cosmic microwave background is predicted by solving the Boltzmann equation for photons in curved spacetime 
(see. e.g., Q). The decoupling of heavy particles from thermal medium in the expanding universe, which is known 
as the generalized Lee- Weinberg problem, is investigated with the use of the Boltzmann equation • The Boltz- 
mann equation for leptogenesis is a recent topics of such early universe physics |J, |5j . Because field theory is more 
0^ ■ fundamental to describe particle processes, it is useful to find the method to derive the Boltzmann equation from the 
field theory. There are two different approaches to derive the Boltzmann equation from the field theory. One is the 
method with the use of the density matrix and the other is the Wigner function approach, which we will adopt in the 
present paper. Thus the Boltzmann equation is regarded as an effective theory of the field theory in a many-particle 
limit, (see e.g. the text book @.) 

On the other hand, recently, the flavor mixing of the neutrinos has been experimentally established by several groups 
0. These results strongly suggest that the neutrinos have non-zero mass and that their interaction basis is not the 
II i mass eigenstate. The propagation of particles with mass mixing might not be trivial in general curved spacetime. 
Actually the gravity effect on the neutrino oscillation has been investigated by several authors ([H, 0, 0, ^| and 
references therein). However, these works focus on how the quantum- mechanical phase of the mixing is affected due 
to the gravity. We might claim that the Boltzmann equation is needed to investigate collective phenomena of many- 
particle systems. The Boltzmann equation of the neutrinos has been investigated by several authors (see [12, Ej, Ej] 
and references therein) . Sigl and Raffelt developed a general formula of the Boltzmann equation for the neutrinos with 
flavor mixing using the density matrix approach |T^ . It has been shown that the same result can be reproduced with 
the Wigner function approach 0, 0] . The Boltzmann equation of the neutrinos is very important to understand 
the early universe. For example, Dolgov et al. have claimed that the flavor equilibrium can be achieved due to the 
flavor mixing effect before the neutrino decoupling, which might lead a significant constraint on the lepton number 
asymmetry of the background neutrinos by combining a constraint from the successful primordial nucleosynthesis 

Era. 

The Boltzmann equation of the neutrinos, in the previous works, have been derived assuming the Minkowski 
background spacetime. The gravity effect on the Boltzmann equation is usually taken into account by replacing the 
Liouville operator in the Minkowski spacetime with that in curved spacetime. The validity of this procedure might 
be worth being checked. The primary aim of the present paper is to derive the Boltzmann equation for fermions 
with flavor mixing in general curved spacetime in a rigorous manner starting from the Dirac equation on curved 
spacetime. For that purpose, we take the covariant Wigner function approach. It is known that a simple Wigner 
function approach encounters the problem of the covariance of the general coordinate transformation. To avoid this 
problem, several authors considered the covariant Wigner function. Winter defined the covariant Wigner function by 
introducing the covariant geodesic distance |17| , while Calzetta, Habib and Hu defined the covariant Wigner function 
using the Riemann normal coordinate |18|. Fonarev developed a very elegant framework for the covariant Wigner 
function using the tangent space. He checked that these three methods yield the locally same equation 0,|2(j. Here 
we adopt the framework by Fonarev to apply it to the fermions with flavor mixing. 

This paper is organized as follows: In section 2 we introduce the covariant Wigner function for mixing fermions, 
and derive its equation of motion. The equation of motion can be regarded as the Boltzmann equation for mixing 
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fermions in curved spacetime. Some details of the derivation are described in the Appendix. In section 3, we show 
that the same equation is reproduced at the lowest order of the expansion with respect to ft in the ultra-relativistic 
limit, which is a useful test of the previous result using the density matrix formalism. In section 4, we demonstrate 
that the familiar formula of the transition probability of the neutrino oscillation is obtained by solving the Boltzmann 
equation. Higher order effects of the ^-expansion on the Boltzmann equation is discussed in section 5. Section 6 is 
devoted to summary and conclusions. 

In the present paper we use the unit in which the light velocity c equals 1. We follow the convention of (H ) 

for the metric g^ v . We use a,P,fj,,y, • • ■ to denote the index of coordinate, r" is the Christoffel symbol, and the 
Riemann tensor is defined R a = d^T^ — d^T^ + r^T^g — r^F^g. In this paper, the superscript, A,B,C ■■■ 
are used to denote the flavor index, while spinor index is omitted. 

II. COVARIANT WIGNER FUNCTION 

In this section we introduce the covariant Wigner function following Fonarev [l^. The covariance for the general 
coordinate transformation is manifest due to the definition with the use of the tangent space pol | . In order to introduce 
the Wigner function, it is instructive to start from considering the energy momentum operator for the fermion fields 



(1) 



where the subscript (/i • • • v) denotes the symmetrization with respect to p and v, A and B denote the flavor index, 
and tr should be taken for the spinor and the flavor indices. The gamma matrix satisfies 



y*7" + 'f'f ' = V"- 
Note that the above expression can be rephrased 



Tuu = -tr 



7( 



(2) 



(3) 



where y„ = y„ — y„. Thus, in general, we consider an operator A which can be written as 

ift^ 



A = -tr 



(4) 



where T consists of 7-matrix and a is a function of —ift\//2. Then we find that A can be written as 



A = -tr 



{nft) 4 



d 4 p^=L= f dS^/^(xje- 2i y a ^/ n a(p u )^ A (x,-y)^ B (x,y), 
V~9W J 



where 



= (l-irV„ + itfVVaV^----)^ A (a:), 
* B {x,y) = (l + y a V a + ^y a y^V a Vp + ---)ij B (x). 



d 4 P 



The expectation value of A is written 

(i)=tr 

J v-sw 

where N AB (x,p) is the covariant Wigner function defined by 

-1 



Ta(p„)N AB (x,p), 



N AB (x,p) = — ' J d A y^f^e- 2i ^ h ^ A {x,-y)^ B {x,y)). 



(5) 

(6) 
(7) 

(8) 
(9) 
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For example the expectation value of the energy momentum tensor is 

(T^) = tr / -fi£= 7ifl p u) N AB (x,p). (10) 

In general curved spacetime, (!),„) diverges. However, how to regularize the divergence is not explicitly specified here 
(see also [U). 

Equation of motion of N AB (x,p) is the (generalized) Boltzmann equation, which we derive from the Dirac equation 
for ij) A (x). In the present paper, we consider the Dirac equation: 

£[^V^ S - M AB (x) - h^J AB {x)]^ B {x) = 0, (11) 

B 

where J AB (x) denotes an effective potential which arises from interaction with background particles |22^, Here we 
adopt the above expression motivated by the interaction of the neutrinos, which are mediated by the vector bosons. 

We briefly summarized the derivation of equation of motion for TV (x,p) in Appendix for being self-contained, 
though the derivation is similar to that in the reference [2£|. Equation of motion for N AB (x,p) is formally obtained 
by the expansion with respect to the power index of Ti. The equation of motion up to the order of Ti 2 is 



h 2 dN AB {x,p) 
16 dp, RvailpG 



VP 



-— V ^R — —N AB (x v) - — R —N AB (x v) ! 

M P H uia ^ ^ ^ (x, p) 4g H a ^ ^ [x, p) ) , 



(12) 



where we defined 



D a = V a + JtuPfiJ- (13) 

OPu 

and a^ p = [7^, 7 p ]/4. We first consider the equation up to the order of % in the next section. The effect of the higher 
order terms in proportion to Ti 2 is discussed in section 5. 



III. EQUATION OF 0{K) 

In this section we focus on the equation of order H: 

J a (p a + y^ a j N(x,p) - \ M{x) + n 1 a J a {x) - ^V a M(x)-^-jN(x,p) = 0, (14) 

where we omit the flavor index of N AB (x,p), for simplicity. Note that we omit the spinor index too. Following the 
procedure in |l3j |. we separate the left handed chirality component from the right handed component. Introducing 
the left- and right- handed chirality projection operators 

1 - 7 5 

Pl = ~^ L , (15) 



1 + T 5 

Pr = —^ L , (16) 



respectively, where 



7 5 = vW^YVyy, (17) 
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we define 

N RL (x,p) = P R N(x,p)P R , (18) 

N L {x,p) = P L N(x,p)P R , (19) 

N R (x,p) = P R N(x,p)P L , (20) 

N LR (x,p) = P L N(x,p)P L . (21) 

In the present paper, in general, we consider written in the form: 

Mx) = PlJl^x) + PrJr„(x), (22) 

where Jl^x) and J R ^{x) are the vector quantities. Projecting Pl and P R from left and right on equation l|14|l . 
respectively, we have 

J a (p a + yA, - hJ Ra {x)\N RL - (m(x) - ^\7 a M(x)-^jN L (x,p) = 0. (23) 
In a similar way, projecting P R and P R from left and right on equation 114f> . we have 

l a (p a + yA* - hJ La {x)\N L - \ M(x) - l ^V a M(x)-^-)N RL (x,p) = 0. (24) 

Combining (|2*5)l and it^l . we eliminate N R l(x,p), and obtain the equation for Nl up to the order of fi 

p a p a - M 2 (x) + ihp a D a + ■ 1 (V a M 2 (x))— - hp a J Lp (x)r^ 

2 dp a 

-hM{x)J Ra {x)Nr\x)pp 1 a 1 p - in{V a M{x))M- 1 {x)p pl a 1 ^N L {x,p) = 0, 

(25) 

where M~ 1 (x) is the inverse of the mass matrix M(x). Similar to the way to derive l|25|l . equations for N R (x,p) can 
be derived, which have the same expression as <|25|) but with replacing R(L) by L(R). 
Now we focus on equation l|25|l . Let us consider a solution of the form: 



N L (x,p)=F^x 1 p)P L Y l . (26) 
Note that Nl(x,p) has the spinor index, but F^(x,p) does not have it. Inserting (|26|l into l|25[) . we finally have 

p a Pa - M 2 (x)+ihp a D a + y(V a M 2 {x))-^F v (x >P ) 

-K a a F v {x,p) +K^F»(x,p) -K v ^{x,p) -ie a ^ v K aP F^{x,p) = 0, (27) 

Kap = hp a J L + hM(x)J Ra (x)M-\x)pi3 + ih(V a M(x))M- 1 (x)p l3 . (28) 



where we defined 



IV. APPLICATION TO NEUTRINOS 

In this section, we apply the result of the previous section to the neutrinos. Setting V Q M (x) = and J Ra (x) = 0, 
equation l|2*7|l reduces to 

(p a Pa - M 2 + ihp a D }jF u - h(p a 3lF v - p a 3 Lv F a + Vv 3 ha F a + ie a ^ vPa J L pF^ = 0. 

(29) 

This is consistent with the result in the reference in which the effect of the gravity is not taken into account. As 
we will show below, in the ultra-relativistic regime, we may set 

F v = p v F, (30) 
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where F — F(x,p) is a scalar function. Substituting l|30[) into (|29|l . we have 

Pv (p a p a - M 2 + ihp a D a - 2hp a J La )F + hp a p a J Lv F = 0. (31) 

The last term of the left hand side of the above equation can be neglected in the ultra-relativistic limit of (p ) 2 3> 
0(M 2 ). In this limit, equation of F is 

(p a p a - M 2 - 2hp a J La + ihp a D a )F = 0. (32) 

Note that F has the flavor index and that F is hermite. Taking the dagger of l|32l) . 

F{p a p a - M 2 - 2hp a J La ) - ihp a D a F = 0. (33) 

Here we assumed j| Q = J La- From (|32[1 and l)33|l we have 

{p a p a - M 2 - 2hp a J La ,F} = 0, (34) 

ihp a D a F = - - [F, M 2 + 2hJ La p a ], (35) 

where {•••,•••} and [•••,•••] denote the anticommutator and commutator, respectively. Equation (|34|l is the constraint 
equation to describe the on-shell condition or the dispersion relation, while H35[l is the Boltzmann equation. This 
result is consistent with the previous result in the flat spacetime limit [oll3lll4| . The gravity effect of this Boltzmann 
equation arises through the differential operator p a D a , with which we can write 



p a D a F{x,p) = p -^J- + Ti vP0 ^-jF(x,p) 



= ±F(x(X),p(X)), (36) 

where A is the afhne parameter. This suggests that the affine parameter will be fundamental to describe the oscillation 
in general curved spacetime. From now on, we use the unit H = 1 for convenience, excepting the case we write h 
explicitly to avoid confusion. 

Here we demonstrate that the familiar formula of the neutrino oscillation probability is obtained by solving the 
Boltzmann equation. For simplicity, we consider the two flavor system with the mass matrix 

M 2 = ( COSd n - Sin / )("$■ ° 2 )( C ° S ° S[n( t ) (37) 

V sin cos J \ m,2 J V — sm 6 cos 1 v ' 

and the source term Ji, a 

•/;:. o 

'La 

In the ultra-relativistic limit, we assume that F might be written 



JLa = ( ) ■ (38) 



f = ( f ;: f ;» ) . (39) 



Writing the off-diagonal component f efl as 

/e M = #e[/ eM ]+i/m[/ eM ], (40) 



(El yields 



where we defined 



p a D a f ee = 2(Im[f efl }, (41) 

p a D a f^ = -2C/m[/ e ^], (42) 

p a D a Im[f etl ] = -C(/ ee - fw) + 2£ite[/ e/ J, (43) 

p a D a Re[f eil ] - -2ZIm[f eiM ], (44) 



i = \(™l- m\) cos20 - \{Jf a - JZ)p a , (45) 
C = j{ml -m?) sin 20. (46) 
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We next demonstrate that the above equations reproduce the well-known formula of the neutrino oscillation on 
the flat spacetime background. We consider one dimensional stationary flux of neutrinos along the z axis, where the 
momentum of the relativistic neutrinos has the component p a = (p, 0, 0, p) . In this case, combining equations Ij41l44|l , 
we have (see also |l3j l 

d 3 1 <9£ d 2 £ 2 + C 2 d 4 d£ Q 2 



and 



dz 



(47) 

{fee + fw) = 0- (48) 



Equation (|48|l yields f ee + =constant, which means the total number of the neutrinos conserves. We must specify 
initial condition to find a solution. We adopt the initial condition so that the electron neutrinos are emitted at z = 0, 
namely, f ee (z = 0) = / ee (0), fnn(z = 0) = 0, and f &il (z = 0) = 0. In the vacuum case J La — 0, integration of l|4Tjl 
yields 



^/ w +^/ w -2^/ ee (0), (40) 



where u> = (m 2 — m\)/2p. Further integration gives the solution 



^-sm^sin 2 ^-^), (50) 



/ee(0) V 4p 

which gives the familiar oscillation probability from the electron neutrino to the mu neutrino. 

The oscillation probability is given by considering the phase evolution of the mass eigenstates in the familiar 
approach. Namely, assuming a wave function that can be decomposed into a coherent superposition of the mass 
eigenstate, the oscillation probability is given by evaluating the transition amplitude from an initial state to other 
state in flavor space. On the other hand, our approach is based on the Boltzmann equation. The equivalence of 
the both approaches is not clearly proved here, however, the above demonstrates the usefulness of the Boltzmann 
approach. Below we brielfy discuss the effect of the gravity on the phase shift by comparing with previous works as 
a test of our framework. 

Consider a trajectory of a particle x a (X) parameterized by the afflne parameter A. Because we are considering the 
ultra-relativistic limit, we assume that the trajectory is a null geodesic. Then, from equation (|36|l . the Boltzmann 
equation in the vacuum is written 

d 3 [ml -m\\ 2 d 

^j(/ee — fufi) + ( ^ J ^X^ ee ~ = ^ 

and 



d_ 

d\ 



{fee + fw) = °' ( 52 ) 



where f ee and / w are understood as the functions along the trajectory f ee — f ee (x(\),p(\)) and = / JliJli (a;(A),p(A)), 
respectively. Along the trajectory we have the formal solution 

/ 2 2 \ 

fee - = Ci + C 2 sin ( ™ 2 ~ mi A(x, p) + 8) (53) 

and 

fee + fw = C 3 , (54) 

where C\, C2, C3 and 8 are the constants which should be determined by initial conditions. Thus the phase of the 
general solution is proportional to the afflne parameter, which may contain the gravitationally induced contribution 
along the geodesic. This result can be consistent with [HJ], in which a resolution for the controversy in the references 
d, 01 is presented (see also 01)- However, our argument assumes that the wave function of the different mass 
eigenstates overlap along the null geodesies. To incorporate the on-shell condition of the different mass eigenstates 
explicitly in our framework, we need to introduce the distribution function using the basis of the mass eigenstate |13| . 
We will revisit this problem beyond the ultra-relativistic limit as a future work. 
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V. HIGHER ORDER EFFECT 

Here wc briefly discuss the higher order effect of the ^-expansion of the Boltzmann equation. As we have shown in 
section 2, the new terms of the order of Ti 2 are 



— 7 /i V«J M (x) — N(p,x) (55) 



T 7-V Q J,(,)- 
and 

h 2 7 a R a [N(p,x)], (56) 

where we defined 

R a [N(p,x)] = l^^R^r^] 

-h?^^ik NM - h^w" M ' (57) 

Similar to section 3, the left- and right- handed components can be separated. Equations corresponding to (|23(l and 
l|2"lf are 

l a (p a + y A.) JVju; + h 2 R a [N RL ] - MN L = 0, (58) 

7 Q (pa + yA« ~ nJ La + ^- V M J ia(^)j|-) n l + ^^R^Nl] - MN RL = 0, (59) 

where we consider the case V a M(x) = and JR a (x) — 0. Combining these equations, 

(p a p a - M 2 + ihp a D a - h 7 a ^p a J Lp ) N L 



+^7 < V (p a (V,J Lp )^--D a {J Lp N L )^ 

-j a 7 PD a D N L + h 2 (j a ^ Pa Rp[N L } + 7 a R a [7%N L ]) = 0. (60) 
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Note that the quantum effect due to spacetime curvature appears as the correction of the order of ft 2 , i.e., expression 
(|56|l . According to the analogy with the equation of the order of h, the effect of the spacetime-curvature alters the 
constraint equation (on-shell condition) [20j , while the variation of the effective potential, 1)55(1 , affects the dynamical 
equation. However, the correction of the higher order terms are small in general cases because the correction is of 
order of Ti/pL, where L is a typical length scale of spacetime curvature or the gradient of the effective potential, and 
p is the momentum of a particle. Thus h/pL is the ratio of the de Broglie wave length to the curvature length scale, 
which is generally small, except for e.g., very early epoch of the universe. 



VI. SUMMARY AND CONCLUSIONS 



In this paper we have derived the quantum kinetic equation for the fermions with flavor mixing in curved spacetime. 
This derivation is based on the covariant Wigner function approach developed by Fonarev |2(j. The result presents 
a rigorous theoretical framework to investigate flavor oscillation phenomena taking the gravitational effect in general 
curved-spacetime into account. The formula is expressed in terms of the expansion with respect to the power index 
of h. The new terms of the order of Ti 2 are the quantum effect due to the gravity, which alter the on-shell condition. 
At present the physical consequence of the correction on the on-shell condition is not clear (cf. see |UI3)- The 
equation of the order of h is consistent with the previous results. We have shown that the equation of the order 
of h reduces to the previous result |12L Il3l with simply replacing the Liouville operator in curved spacetime 
to that in the Minkowski spacetime. It is also demonstrated that the familiar formula for the vacuum neutrino 
oscillation probability can be obtained by solving the Boltzmann equation. Then it is shown that our approach gives 
the consistent results which relies on the phase evolution. However, these results assume the ultra-relativistic limit. 
We have also derived the general Boltzmann equation which does not assume the ultra-relativistic limit in section 
3. Analysis of the equation will clarify the effects of the finite mass on the neutrino's propagation in general curved 
spacetime. We will revisit this problem in future work. 
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APPENDIX A: MATHEMATICAL FORMULAS 



Here we derive equation for the covariant Wigner function N AB defined in The derivation is similar to that in 
the reference |'2('| . excepting our generalization of the mass term and the flavor index, however, we write this Appendix 
for being self-contained. It is useful to define the derivative operator: 



v Q = v Q -r^^, (Ai) 



d_ 

dp' 

then we find 

V Q / = 0. (A2) 

Equations JJjJ and JJJ) are, respectively, written as 

* A (x, -y) = exp(~y a \7 a )^ A (x) (A3) 

and 

* B (x,y)=exp(y a V a )ij B (x), (A4) 
For arbitrary operators A and the identities 

oo 1 

[A,e B ] = -J2-[B,[---AB,A]--.}}e B (A5) 

n=l 

and 

v n— 1 v 7 7 

are satisfied. Using these relation we have 

V a * A (x,±y) = e^-V^z) - H a (x, ±y)* A (x, ±y) (A7) 

and 

±y) = ±V a V A (x, ±y) ± G a (x, ±y)V A (x, ±y), (A8) 



where we defined 



and 



H a (x, ±y) = —TV" 1 ■■■V Vn ^ ,[•", [V„ B , V a ] ■ ■ •]] (A9) 



00 \n 

G a (x, ±y) = Yl T^TlV^ 1 '" yUn [ ^ U1 ' [ " ' ' [Wv " ^ a] " ' ]] ' (A10) 

n— 1 ^ 



9 



respectively. We also have 

e ±v"v„J2(M AB (x) + h^J AB (x))i> B (x) 



= J2(M AB (x) + h^J AB (x)) * B (x, ±y) + J2 L AB {x, ±y)^ B (x, ±y), (All) 

B B 

with defining 

L AB (x, ±y) = £ ^pV 1 ■ ■ ■ V^u, ■ ■ • V„„ (m ab (x) + h^J AB {x)) . (A12) 



n=l 



Next we consider D a N AB (x,p), where D a is defined by l|13|) . Using eas. l|A7|) and (|A8|I . we have the expression: 
D a N AB (x,p) = 2 i Pa N AB {x,p)-^f^ (d'ye^y^ 



n*"" ^'^ inn) 4 , 

x (2(e-^ V a <p A (x))* B (x, y) - G AB {x, y)), (A13) 

where 

Gf^y) = V A (x,-y)(G a (x,y)y B (x,y)) - (G a (x,-y)y A (x,~y))y B (x,y) 

+2(H a ( Xl -y)^ A {x 7 -y)) y). (AH) 



Projecting ih r y a (x)/2 on the both side of (|A13|I . we have 



xr(x)^(e-y'^^ a ^ A (x))^ B (x,y)-G AB (x,y)). (A15) 

Using the Dirac equation pi(l and 

V^O) = 0, (A16) 



we have 



l a (x) (p* + ^D a )N AB (x,p) =Y,M AC (x, J±^L)N CB (x,p) 



\f-g{x) 

' 2 (irk) 4 



d*ye- 2i y *>e/ h r(x)(G AB (x,y)), (A17) 



where 



M AB [x,^^-) = Al AB (x)+h^J AB 



(A18) 



Then we focus on the last term of the right hand side of (|A17|) . Using to denote the binomial coefficient, one can 
derive the following recursion formula: 

G A{n) = -A^- 1] ^ A {x, -y) + i^^tf A (x, -y) - £ C^ 2 R^~ k > G A{k \ (A19) 

^ fc=i 



10 



where we defined 



G Mn+i) = yVl . . . yVn ^ > . . s [v^, V Q ] • • -y) (n > 1), 

G£W=V a * A (:r,-y), 



(A20) 



(«>3), 



(A21) 



with 



A( n ) — „«1 . . . q ,Vn A 
4(1) — v "iA 



(n>2), 



(A22) 



which satisfies the commutator relation [V Q) Vp]i() A (x) = A a pip (x). In a similar way, we also have the following 
recursion formula 



with 



GfW = 9 B (x,y)A<£-» + R^>^ B (x,y) - £ C^ 2 R.£~ k > G B W , 

y k=i 



G B(n+1) = y^--- y v - [V^ ,[..-, [V„ B , VJ • • -p B (.T, y) (n > 1) 
G B U=V a * B (x,y). 



(A23) 



(A24) 



A Repeated use of the recursion formula (|A19I) or l|A23|) allows one to rewrite G a (x, — y)$? A (x, — y) and 
H a (x, —y)^> A (x, —y) or G a (x, y)^> B (x, y) in equation (|A14() . in terms of d , i> A (x, — y)jdy v , (x, — y) and 

^ A (x, —y), eliminating the derivative terms higher than the second derivatives. Then the formal expression for 
G AB is written in the form, 



G AB (x,y) 



y K^Rf v v v y K^RfvJL 

' < (k)a p " / j (k)a P Qy V 

k=2 k=2 y 

oo 

-yK[lf a Af^^ A (x,^ B (x,y)) 

k=2 

oo 

^ A (x,-y)* B (x,y))yKlXA { r\ 



(* A (x,-y)<i> B (x,y) 



k=2 



where 



K (i)P _ xP C (i) i V V • • • V G (l) #(*»)&> . . . ??( fcl 



)/3 



(A25) 



(A26) 



n=l fci=2 fc„=2 



and C$ and C(^ fcl ... fc are numerical coefficients, which can explicitly be found using generating functional method 
[20|. We omit the general expression, however, the lowest order coefficients are, for example, 



(1) _ 




-1)* 


(fc) _ 


k\ (k + 


1)! 


(2) _ 


(-1)* l + (- 


-l) k 


(fc) 


fc! (fc + 


1)! 


(3) 


(-l) fc (2fc-l) 




(fc) ~ 


kl 




(4) 


1 




(fc) ~ 


fc!' 





(A27) 

(A28) 

(A29) 
(A30) 
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Using l|A25(l . we have 



where 



jh V-g( x ) 

2 {wh) 4 

oo 

E 

k=2 

ih , 



ih . 



-(4)/3 .(fc-1) 



r «/3 _ rf?^) i V V ... V <7 W T?(*n)i8n ... 7? 



(fei)/3 
'/9a ' 



n=lfci=2 fc„=2 



2 / * , fiai^;i/3'"i/js 



" Vk dp Vl dp v J 



ih 



A 



(A31) 

(A32) 

(A33) 
(A34) 



Finally we find that equation N (x,p) follows 



ih 



7 a (*) Pa + ~7j~D a )N AB (x,p) -Y^M 



ACf n 



-ih d 
~~dp~ 



N CB (x,p) 



'(*)£ 



-(2)/3 1? (fe) I / * ?i r (l)/5 7? (fc)f n i v.iB 

-(jfe)c« /9 ^ + 2 0)« (3 

fe=2 L x 



-(4)0 „(fe-l) 



(A35) 



[1] W. Hu and S. Dodelson, Annual Review of Astronomy and Astrophysics, 40, 171 (2002) 
[2] J. Bernstein, Kinetic Theory in the Expanding Universe, (Cambridge University Press, 1988) 
[3] E. W. Kolb and M. S. Turner, The Early Universe, (Addison- Weesley Publishing, 1990) 
[4] W. Buchmuller and M. Plumacher, Int. J. Mod. Phys. A15, 5047 (2000) 

[5] T. Endoh, S. Kaneko, S. K. Kang, T. Morozumi, M. Tanimoto, Phys. Rev. Lett. 89, 231601 (2002) 

[6] S. R. de Groot, W. A. van Leeuwen, and Ch. G. van Weert, Relativistic Kinetic Theory, (North-Holland Publishing, 1980) 

[7] e.g., S. Fukuda, et al. Phys. Rev. Lett. 85, 3999 (2000); Q. R. Ahmad, et al. Phys. Rev. Lett. 89, 1301 (2002) 

[8] D. V. Ahluwalia and C. Burgard, Gen. Rel. Grav. 28, 1161 (1996); D. V. Ahluwalia and C. Burgard, Phys. Rev. D 57, 

4724 (1998); Y. Kojima, Int. J. Mod. Phys. Lett. A 11 2965 (1996) 
[9] T. Bhattacharya, S. Habib, and E. Mottala, gr-qc/9605074 T. Bhattacharya, S. Habib, and E. Mottala, Phys. Rev. D 59, 
067301 (1999) 

[10] J. Wudka, Phys. Rev. D 64, 065009 (2001) 

[11] K. Konno and M. Kasai, Prog. Theor. Phys. 100, 1145 (1998) 

[12] G. Sigl and G. Raffelt, Nucl. Phys. B406, 423 (1993) 

[13] M. Sirera and A. Perez, Phys. Rev. D59, 125011 (1999) 

[14] S. Yamada, Phys. Rev. D62, 093026 (2000) 

[15] A. D. Dolgov, S. H. Hansen, S. Pastor, S. T. Petcov, G. G. Raffelt, D. V. Semikoz, Nucl. Phys. B 632, 363 (2002) 
[16] A. D. Dolgov, Physics Report, 370, 333 (2002) 
[17] J. Winter, Phys. Rev. D 32, 1871 (1985) 

[18] E. A. Calzetta, S. Habib and B. L. Hu, Phys. Rev. D 37, 2901 (1988) 
[19] O. A. Fonarev, Soviet Phys. J. 33, 759 (1990) 
[20] O. A. Fonarev, J. Math. Phys. 35, 2105 (1994) 

[21] N. D. Birrell and P. C. W. Davies, Quantum fields in curved space, (Cambridge University Press, 1982) 

[22] D. Notzold and G. Raffelt, Nucl. Phys. B307, 924 (1988) 

[23] O. A. Fonarev, Phys. Lett. 190A, 29 (1994) 

[24] K. Pirk, and G. Boerner, Class. Quant. Grav., 6, 1855 (1989) 



